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We propose a framework to model elastic properties of polycrystals by coupling crystal orientational degrees 
of freedom with elastic strains. Our model encodes crystal symmetries and takes into account explicitly the strain 
compatibility induced long-range interaction between grains. The coupling of crystal orientation and elastic 
interactions allows for the rotation of individual grains by an external load. We apply the model to simulate 
uniaxial tensile loading of a 2D poly crystal within linear elasticity and a system with elastic anharmonicities 
that describe structural phase transformations. We investigate the constitutive response of the poly crystal and 
compare it to that of single crystals with crystallographic orientations that form the polycrystal. 
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A study of the mechanical properties of polycrystals is im- 
portant as most technologically important materials exist in a 
polycrystalline state. A polycrystal is an aggregate of grains 
that have different crystallographic orientations. The prop- 
erties of a polycrystal depend on its texture (distribution of 
crystallographic orientations). It is desirable to understand 
how the strains are distributed due to an applied external load 
and how this influences the average elastic moduli of the 
polycrystaltl. An important feature of some polycrystalline 
metallic alloys and ceramics is the presence of domain walls 
within the grains due to an underlying martensitic transition. 
This microstructure influences the response of the material to 
an external load. For example, in martensites, the deformation 
is accompanied by the motion of the domain walls. Texture 
evolution caused by rotation of grains under the application of 
an external load is another important factor that influences the 
mechanical behavior of a polycrystal. Polycrystalline spec- 
imens exhibit significant grain rotation in the plastic regime 
to accommodate crystallographic slipB. Grain rotations to 
~ 1° have also been observed in atomistic simulations™ of 
nanocrystals, even at low strains of ~ 4%. 

The problem of finding the effective properties of polycrys- 
tals has been studied by analytical methodsElB. However, the 
complex geometry of a polycrystal or the long-range elastic 
interactions between the grains are often not accounted for in 
these methods. In fact, these approaches are only able to give 
bounds on the effective properties. Recently, the mechani- 
cal deformation -of polycrystals has been studied by atom- 
istic simulationsm which have been limited to nano- sized 
grains. Simulating bulk systems with atomistic simulations 
requires enormous computational power and hence contin- 
uum simulations that can cover a range of intermediate length 
scales are essential for describing the microstructure at the 
sub-micron scale. Several phase-field models have been pro- 
posed to model grain growth phenomenaBM Although these 
models correctly describe the grain morphologies and the do- 
main growth laws, the issues of elasticity and material specific 
crystal symmetries are usually not addressed. Recently Elder 
et alH studied elastic and plastic effects using a model for- 
mulated to describe pattern selection. The model, in terms 
of particle density fields, is specific only to certain symme- 
tries that are selected by appropriate choice of wavelengths 
in the free energy, which does not contain experimentally 



measured quantities such as elastic constants. In this Letter, 
we propose a polycrystal model based on continuum elastic- 
ity that can be applied to any crystal symmetry and has the 
appropriate single-crystal elastic constants as input parame- 
ters. In this model, elastic strains are coupled to a phase field 
model through an orientation field that is determined from a 
multi-component order parameter describing the crystal ori- 
entations. Due to this coupling, the strains in each grain as 
well as the grain orientations can change under an external 
load. This experimentally relevant feature is not accounted 
for in models, that consider static grains created by Voronoi 
constructionH In the present work, we determine the me- 
chanical properties of linear elastic materials and those de- 
scribed by nonlinear elasticity, such as martensites. 

The free-energy functional is written as F = F grain + 
F elastic + Fioad, where Fg ra i n is the free energy density due 
to the orientational degrees of freedom of the polycrystal, 
Feiastic represents the elastic free energy and Fi oa d is the free 
energy contribution due to an external applied load. The poly- 
crystalline system is described by a set of Q non-conserved 
order parameters^ (771, 772, Vq)- A given grain orientation 
corresponds to one of the Q order parameters being positive 
nonzero while the rest are zero. The free energy F grain is 
given by 
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For ai, a2 < and 0,3, > 0, the first two terms in equa- 
tion (1) describe a potential with Q degenerate minima corre- 
sponding to Q grain orientations. The gradient energy (K > 
0) represents the energy cost of creating a grain boundary. It is 
also possible to associate an orientational field 6(ff, r), where 
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Thus, there are Q orientations between and m . For the 
elastic free energy the linearized strain tensor in a global 
reference frame is defined by = (uij + Uj^)/2 (i = 
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1,2: j = 1,2), where Ui represents zth component of the 
displacement vector and Uij is its jth displacement gradi- 
ent. For illustration, we consider a 2D lattice with square 
symmetry and use the symmetryj^dapted linear combinations 
of the strain tensor defined astil e\ = (e xx + e yy )/V2, 



damped dynamics for the strains is 
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To generalize this where 72 and 73 are the appropriate dissipation coefficients 



theory for the case of a polycrystal, the strain tensor in 
a rotated frame is calculated as R(0(ff))eR T (9(ff)), where 
R{Q(rjj) is a rotation matrix. Using this transformation, 
the elastic free energy in a global frame of reference is 
Feiastic = J dr{^ ei 2 + ^e 2 2 + 4^e 3 2 + f n i(e u e 2 , e 3 ) + 
^■(Ve2) 2 + ^-(Ves) 2 } where ei, e 2 , e 3 are defined as 
ei = ei, e 2 = e 2 cos[20(t7)] + >/2e 3 sin [20 (77)] and e 3 = 
-(l/y/2)e 2 sin [20 (77)] + e 3 cos [20 (77)]. The orientation field 
0(77) is determined from the minima of free energy in (1) using 
(2). Here A 1 = C n + C 12 , A 2 = C n - C 12 and A 3 = 4C 44 , 
where Cn, C12 and C44 are the elastic constants for a crys- 
tal with square symmetry. K 2 and K3 are the appropriate 
gradient coefficients that in principle can be obtained from 
experimentally measured phonon dispersion data. The term 
fni(ei,e2, e 3 ) represents the nonlinear part of the elastic free 
energy and is crucial in describing structural phase transitions. 

In this work, we are interested in simulating a uniaxial load- 
ing experiment. If we choose the x axis to be the loading di- 
rection, the free energy contribution due to the external load is 
Fioad = - Jdfae xx = - J dr^(e 1 +e 2 ) = - f dr-^(e 1 + 

e 2 cos [20 (77)] — V^es sin [20(77)]). The strains e\ , e 2 and^ 
are not independent but satisfy a compatibility relationshipt3: 
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^"a?^-v"e; e 3 = - Usin S a method in- 
duced earlier for single crystal martensitic transformationsuJ, 
the strain e\ may be eliminated using compatibility, to express 
the effective free energy F eff = F e i astic + Fi oad as 
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where T 2 (k), T 3 (/c) represent Fourier transforms of 
e 2 cos[20(77)] - V^e3sin[20(77)] and e 2 (sin [20(77)] /V^) + 
e 3 cos [2 0(77)] respectively, C2(fc) = (k x 2 — k y 2 ) / (k x 2 -\- k y 2 ) 
and C 3 (fc) = VSk x k y / (k x 2 + /c^ 2 ). The long-range terms 
ensure that compatibility is satisfied within the grains as well 
as at the grain boundaries. 

The dynamics of the grains is given by Q equations 
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where 7^ is a dissipation coefficient and i = 1, cor- 
respond to Q grain orientations. The corresponding over- 



for the strains and F = F gra i n + F e / / is the total free energy 
of the system. 

The mechanical properties of many materials are well de- 
scribed by the harmonic approximation for which the non- 
linear term f n i(ei, e 2: e 3 ) = 0. For a homogeneous single 
crystal e 2 , e 3 and 77 are constant and for Cu A\ = 289.8 
GPa, A 2 = 47.0 GPa and A 3 = 301.6 GPa ( Cn = 168.4 
GPa, C12 = 121.4 GPa, C 44 = 75.4 GPa). For the pa- 
rameters in F grain we choose a\ — a 2 = —A 2 , a 3 = A 2 , 
C14 = 2A 2 , Q = 5 and m = 45°. We choose the gradi- 
ent coefficients in terms of an arbitrary length scale 5 so that 
K = K 2 = K 3 = A 2 S 2 and lengths are scaled by r = 5£. 
The free energy in (1) then has five degenerate minima defined 
by o (ff) = 0°, 11.25°, 22.5°, 33.75°, 45°, corresponding to 
five different grain orientations. 

To study the polycrystal, we first generate an initial poly- 
crystalline configuration by solving (4) and (5) with a = 0, 
using random initial conditions. For all simulations in this pa- 
per, we assume 7^ = 71 = 72 = 7 and use rescaled time 
t* = t \A 2 1 7. For a = 0, the elastic effects do not influ- 
ence the grain growth as all the strains vanish. Grains with 
orientations O (77) = 0°, 11.25°, 22.5°, 33.75°, 45° nucleate 
and coarsen. We arrest the system in a given configuration by 
suddenly changing the value of the parameter a\ from — A 2 
to — 16 A 2 . This increases the free-energy barriers between 
the crystalline states and the growth stops. With the arrested 
polycrystal configuration as the initial condition, we simulate 
a quasi- static uniaxial tensile loading using (4) and (5). The 
stress a is varied in steps of 0.06 GPa and we let the strains 
relax after each change for t* = 25 steps. Figure 1(a) shows 
the spatial distribution of the polycrystal orientation a (f) at 
a loading of a = 2.35 GPa for a system of size 1285 x 1285. 
We note that there is no significant motion of the grain bound- 
aries from the initial arrested configuration to the configura- 
tion depicted in Fig. 1(a) (the individual grains have rotated 
by a small amount (~ 0.01°), consistent with the coupling 
between stress and the orientation). In Fig. 1(b), we show the 
corresponding distribution of uniaxial strain e xx . The strain 
distribution is anisotropic as e xx in a grain depends on the ori- 
entation. 

In Fig. 2, we show the variation of the average strain (e^) 
with the load a. For comparison, we also plot the analo- 
gous single crystal curves with crystallographic orientations 
that constitute the polycrystal configuration in Fig. 1 (single 
crystal simulations were performed using only one orientation 
but with identical free energy parameters and loading rate as 
the polycrystal). The Young's modulus of the simulated poly- 
crystal was ~ 126 GPa. This is in. the range of experimen- 
tally measured values of 124 GPaH and 129.8 GPaLJ quoted 
for bulk polycrystalline Cu. The result is not sensitive to the 
choice of parameters for the polycrystal phase field model, at 
least in the linear elastic regime. 



3 



Another important class of materials that can be studied us- 
ing this approach are martensites that undergo a displacive 
structural phase transformation. The transformed phase is 
characterized by a complex arrangement of crystallographic 
variants known as twins. We consider the case of a 2D square 
to rectangle transition for which the deviatoric strain 62 is the 
appropriate order parameter. For the high temperature square 
phase e2 = and for the low temperature martensitic phase 
e2 = ±eo, corresponding to the two rectangular variants. 
This system exhibits the so called shape-memory effect which 
is governed by the motion of martensitic domains. The mi- 
crostructure depends on the underlying crystal symmetry and 
hence the displacements and domain wall orientations of the 
atoms in each grain depend on the grain orientation. Thus 
the shape memory effect will be influenced by the texture and 
hence it is important to compare the mechanical response of 
single and poly crystal martensites. 

The anharmonic contribution to the elastic free energy is 
given by f n \ = f e2 4 + f e2 6 and describes a first order tran- 
sition for a < 0. We choose A 1 = 140 GPa, A 3 = 280 
GPa, a = -1.7 x 10 4 GPa^id p = 3 x 10 7 GPa, parame- 
ters that correspond to FePdlij. The constant A2 depends on 
the temperature and we choose A2 = — 3 GPa, a tempera- 
ture in the martensitic phase. The measured gradient coeffi- 
cient K2/CL0 2 = 25 GPa, where ao is the lattice spacing of 
the crystal and assume K3 = since the deviatoric strain 
is the dominant mode of deformation. The parameters for 
Fgrain are ai = a 2 = -\A 2 \, a 3 = \A 2 \, a 4 = 2\A 2 \ and 
the grain boundary coefficient is chosen in terms of the lat- 
tice spacing ao to be K/clq 2 = 10 4 | A2 1 . The space variable 
is rescaled by introducing a dimensionless length scale ( so 
that r = (100ao)C The maximum orientation is chosen to 
be m = 30° so for Q = 5, the allowed orientations are 0°, 
7.5°, 15°, 22.5° and 30°. Employing the above set of param- 
eters, we simulate the martensitic domain structures using (4) 
and (5). With same procedure as for the linear elastic case, 
a stable polycrystal configuration is obtained. When the ap- 
plied stress a = 0, the parameter A2 = — 3 GPa ensures that 
the system is well in the martensite phase and domains of the 
two rectangular variants (twins) are formed. After obtaining 
a stable martensitic polycrystal, the loading process is simu- 
lated by quasi- statically varying the stress in steps of 0.38 GPa 
upto a maximum stress 15 GPa (after each stress change, the 
system is allowed to relax for t* = 25 steps). The system is 
then unloaded by decreasing the stress to zero at the same rate. 

Figure 3 shows the evolution of variants and the grains at 
different stress levels during the loading-unloading process 
for a system of size 12800ao x 12800ao. The left column 
shows the distribution of e 2 (r) (deviatoric strain in a global 
frame) and the right column shows the corresponding distri- 
bution of the orientations 0(f). It is clear from Fig. 3 that the 
domain wall orientations depend on the orientation of a grain. 



On loading, the simulated polycrystal starts to detwin (favored 
variants grow at the expense of unfavored ones). However, 
even at the maximum load of a = 15 GPa, some unfavorable 
variants persist. On complete unloading, a domain structure 
is nucleated again due to inhomogeneities in the polycrystal. 
However, this domain structure is not the same as that before 
the loading indicating an underlying hysteresis. The orien- 
tation distribution is also influenced by the external load, as 
can be seen in the right column of Fig. 3. The grains with 
large misorientation with the loading axis rotate significantly 
(~ 10°) while the grains having lower orientation do not ro- 
tate as much. The mechanism of this rotation is the desire of 
the system to maximize strain in the loading direction so as to 
minimize the elastic free energy. At high stress, some grain 
boundaries start moving to accommodate the applied stress, 
as is clear from the orientation distribution at a = 15 GPa. 

The stress-strain curve corresponding to Fig. 3 is shown 
in Fig. 4. Also shown are single crystal curves for all five 
orientations that constitute the polycrystal of Fig. 3. We 
observe that the hysteresis for the polycrystal case is much 
smaller than that for a single crystal oriented along the loading 
axis. Our findings regarding the hysteresis are consistent with 
the fact that polycrystals have poor shape memory properties 
compared to single crystals!! The simulations also indicate 
that grain rotations will influence the mechanical properties 
of shape memory alloys. Recently, in situ measurements of 
texture evolution during compression experiments on Ni-Ti 
shape memory alloystZl have been reported. However, these 
experiments cannot predict whether the changes in texture are 
due to detwinning or rotation of grains. Our simulations show 
that both these processes can contribute to texture evolution. 

In summary, we have proposed a framework to study the 
mechanical properties of polycrystals in which the long-range 
elastic interaction between grains and the connectivity of the 
microstructure is taken into account. The approach can be 
extended to any crystal symmetry or loading (e.g., shear) and 
does not require any a priori assumption of grain shapes or the 
microstructure. An important feature of our work is the cou- 
pling between the grain orientation and elasticity. We have ap- 
plied the model to study mechanical properties of linear elastic 
and martensitic materials. For the linear elastic case the ob- 
served grain rotations are small (~ 0.01°) and hence do not 
influence the mechanical properties. In contrast, the marten- 
sitic case shows significant grain rotation (~ 10°) due to ac- 
commodation of the transformation strain. This behavior is 
sensitive to the choice of parameters of the polycrystal model 
(energy barriers between grains) and therefore determination 
of these parameters from experiment or atomistic simulations 
will allow accurate prediction of mechanical properties. 
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Figure Captions: 

Figure 1: Spatial distribution of orientation angle 0(f) 
(snapshot (a)) and uniaxial strain e xx (f) (snapshot (b)) for 
stress a = 2.35 GPa. 

Figure 2: Average uniaxial strain (e xx ) for the linear elastic 
case as a function of the load a. The curves correspond to a 
polycrystal (o) and single crystals with Oq — 0° (x), Oq — 
11.25° (+), 6 = 22.5°(*), 6 = 33.75°(D) and <9 = 45° 
(o). 

Figure 3: Spatial distribution of the deviatoric strain in a 
global frame, 62 (r), (snapshots (a),(c),(e) and (g)) and orien- 
tation angle 9(f) (snapshots (b),(d),(f) and (h)). The corre- 
sponding stress levels are a = ((a) and (b)), a = 2.69 GPa 
((c) and (d)), a = 15 GPa ((e) and (f)) and a = (after un- 
loading) ((g) and (h)). 

Figure 4: Average uniaxial strain (e^) for the martensite as 
a function of the load a. The curves correspond to a polycrys- 
tal (o) and single crystals with = 0° (x), = 7.5°(+), 
6> = 15°(*), 6 = 22.5°(D) and <9 = 30°(o). 
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